subroutine ccsdt(icore, eccsdt) 
use mod_iop
use mod_orbit
use mod_symm
use mod_ioff
use mod_size 
use mod_t
implicit none 
integer prec2, isd2, igpu2, totgpu2
common/precsoc/prec2, isd2, igpu2, totgpu2
!
real*8, intent(in) :: eccsdt
!
real*4, intent(in) :: icore(*)
!
real*4, allocatable,device:: icorg(:)
!
integer ioffabci1(noa, 8)
real*8 e4t, e5st, tt(10), rt(10), eccsd, escf
integer it1, it1aa, it1bb, it2, it22, iijka, iabij, iend, iabci
integer itmp1, itmp2, itmp3, itmp4, itmp5, itmp6, itmp7
integer irrep, ioffie, ie, icrsizg
integer it21, it23, it221, it222, it223, iijka1, iijka2, iijka3
integer iabij1, iabij2, iabij3, iabci1, iabci2, iabci3 
integer nsize
integer, external :: irpdso, idsymsoc
!
if(totgpu==1) then 
   icrsizg = igpusize!10GB
   write(6,*)'memory allocated in GPU', icrsizg/1024/1024*4, 'MB'
   allocate(icorg(icrsizg))
endif
!
if(prec==8) then 
  call load(icore, icore)
elseif(prec==4) then 
  call loads(icore, icore)
endif
!
do irrep = 1, nirrep
   nvvaa(irrep) = irpdso(irrep, vrta, vrta, 0)
   nvvbb(irrep) = irpdso(irrep, vrtb, vrtb, 0)
   nvvab(irrep) = irpdso(irrep, vrta, vrtb, 1)
   nvvaa2(irrep) = irpdso(irrep, vrta, vrta, 1)
   nvvbb2(irrep) = irpdso(irrep, vrtb, vrtb, 1)
!
   nooaa(irrep) = irpdso(irrep, popa, popa, 0)
   noobb(irrep) = irpdso(irrep, popb, popb, 0)
   nooab(irrep) = irpdso(irrep, popa, popb, 1)
!
   nvoaa(irrep) = irpdso(irrep, vrta, popa, 1)
   nvoab(irrep) = irpdso(irrep, vrta, popb, 1)
   nvoba(irrep) = irpdso(irrep, vrtb, popa, 1)
   nvobb(irrep) = irpdso(irrep, vrtb, popb, 1)
enddo
ie = 1 
it1 = ie + (nva + noa + nvb + nob)*isd
it2 = it1 + nifvo*isd 
it22 = it2 + nvvoo*isd
iijka = it22 + nvvoo*isd
iabij = iijka + nvooo*isd
itmp1 = iabij + nvvoo*isd
itmp2 = itmp1 + nvvv*isd
itmp3 = itmp2 + nvvv*isd 
itmp4 = itmp3 + nvvv*isd 
itmp5 = itmp4 + nvvv*isd 
iend  = itmp5 + nvvv*isd
!
iabci = iend 
nsize = max(idsymsoc(1, vrta, vrta, 0, vrta, popa, 1), & 
            idsymsoc(1, vrtb, vrtb, 0, vrtb, popb, 1))
iend = iabci + nvvvo*isd
!write(6,*)'tmp', isd, nvvoo, nvvv, nifvo
write(6,*)'memroy needed for ccsdt aaa and bbb',iend/1024/1024*4, 'MB'
e4t = 0.d0; e5st = 0.d0 
e4ttest = 0.d0; e5sttest = 0.d0 
tgen = 0.d0; tvpq = 0.d0; ttmpe = 0.d0; ttmpe2 = 0.d0
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=============for aaa calculation============!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
call cpu_time(tt(1)); call walltimeso(rt(1))
!ldebug=.false.
!ldebug=.true.
if(prec==8.and.totgpu==0) then 
   call ccsdt_sym_aaa_d(icore(ie), icore(it1), icore(it2), &
                   icore(it22), icore(iijka), &
                   icore(iabij), icore(itmp1), icore(itmp2), & 
                   icore(itmp3), icore(itmp4), &
                   icore(iabci), icore(iend), e4t, e5st)   
elseif(prec==4.and.totgpu==0) then 
   call ccsdt_sym_aaa_s(icore(ie), icore(it1), icore(it2), & 
                   icore(it22), icore(iijka), &
                   icore(iabij), icore(itmp1), icore(itmp2), & 
                   icore(itmp3), icore(itmp4), &
                   icore(iabci), icore(iend), e4t, e5st)   
elseif(prec==8.and.totgpu==1) then
   write(6,*)'================calculate aaa in gpu. dp================='
   if(iend > icrsizg) then 
       write(6,*)' attention memory not enough for abciaaa in gpu'
       write(6,*)'integrals would be divided into small part'
       call ccsdt_sym_aaa_dg_l2(icore(ie), icore(it1), icore(it2), &
                   icore(it22), icore(iijka), &
                   icore(iabij), icore(itmp1), icore(itmp2), & 
                   icore(itmp3), icore(itmp4), &
                   icore(iabci), icore(iend), e4t, e5st, icorg)   
  else 
  igpu=1; igpu2=1
  call myicopy_g(icore(ie), icorg(ie), nva+noa+nvb+nob)
  call ccsdt_sym_aaa_dg(icorg(ie), icorg(it1), icorg(it2), &
                  icorg(it22),icorg(iijka), &
                  icorg(iabij), icorg(itmp1), icorg(itmp2), &
                  icorg(itmp3), icorg(itmp4), &
                  icorg(iabci), icorg(iend), e4t, e5st, icore(it1))
  endif
elseif(prec==4.and.totgpu==1) then
   write(6,*)'================calculate aaa in gpu. sp================='
  if(iend > icrsizg) then
       write(6,*)'attention!!!! memory not enough for (T) calculations'
       write(6,*)'integrals would be divided into small part'
       call ccsdt_sym_aaa_sg_l2(icore(ie), icore(it1), icore(it2), &
                   icore(it22), icore(iijka), &
                   icore(iabij), icore(itmp1), icore(itmp2), & 
                   icore(itmp3), icore(itmp4), &
                   icore(iabci), icore(iend), e4t, e5st, icorg)   
  else 
  igpu=1; igpu2=1
   call myicopy_g(icore(ie), icorg(ie), nva+noa+nvb+nob)
  call ccsdt_sym_aaa_sg(icorg(ie), icorg(it1), icorg(it2), &
                  icorg(it22),icorg(iijka), &
                  icorg(iabij), icorg(itmp1), icorg(itmp2), &
                  icorg(itmp3), icorg(itmp4), &
                  icorg(iabci), icorg(iend), e4t, e5st, icore(it1))
  endif
endif
write(6,*)'e4t and e5st test for aaa', e4ttest, e5sttest
write(6,'(A, F8.2, F8.2, F8.2)')'tgen,tvpq and ttmpe', tgen, tvpq, ttmpe
write(6,'(A, F8.2)')            'ttmpe2             ', ttmpe2
call cpu_time(tt(2)); call walltimeso(rt(2))
if(prec==8) then 
  call load(icore, icore)
elseif(prec==4) then 
  call loads(icore, icore)
endif
ie = 1 
it1aa = ie + (nva + noa + nvb + nob)*isd
it1bb = it1aa + irpdso(1, vrta, popa, 1)*isd
it21   = it1bb + irpdso(1, vrta, popa, 1)*isd
it22   = it21   + nvvoo*isd
it23   = it22   + nvvoo*isd
it221  = it22   + nvvoo*isd
it222  = it221  + nvvoo*isd
it223  = it222  + nvvoo*isd
iijka1 = it222  + nvvoo*isd
iijka2 = iijka1 + nvooo*isd
iijka3 = iijka2 + nvooo*isd
iabij1 = iijka2 + nvooo*isd
iabij2 = iabij1 + nvvoo*isd
iabij3 = iabij2 + nvvoo*isd
itmp1  = iabij2 + nvvoo*isd
itmp2  = itmp1  + nvvv *isd
itmp3  = itmp2  + nvvv *isd
itmp4  = itmp3  + nvvv *isd
itmp5  = itmp4  + nvvv *isd
itmp6  = itmp5  + nvvv*isd
iend   = itmp6  + nvvv*isd
!
iabci1 = iend 
iabci2 = iabci1 + nvvvo*isd
iabci2 = iabci1 + idsymsoc(1, vrta, vrta, 0, vrta, popa, 1)*isd
iabci3 = iabci2 + nvvvo*isd
iend   = iabci2 + nvvvo*isd
!
write(6,*)'memroy needed for ccsdt: aab and bba ', iend/1024/1024*4, 'mb'
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=============for aab calculation============!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
tgen = 0.d0; tvpq = 0.d0; ttmpe = 0.d0; ttmpe2 = 0.d0
call cpu_time(tt(3)); call walltimeso(rt(3))
!ldebug=.true.
if(prec==8.and.totgpu==0) then 
  call ccsdt_sym_aab_d(icore(ie), icore(it1aa), &
            icore(it21), icore(it22),  &
            icore(it221), icore(it222),  &
            icore(iijka1), icore(iijka2), &
            icore(iabij1), icore(iabij2),  &
            icore(itmp1), icore(itmp2), icore(itmp3), &
            icore(itmp4), icore(itmp5), icore(itmp6), &
            icore(iabci1), icore(iabci2), &
            icore(iend), e4t, e5st)
elseif(prec==4.and.totgpu==0) then 
  call ccsdt_sym_aab_s(icore(ie), icore(it1aa),  &
            icore(it21), icore(it22),  &
            icore(it221), icore(it222),  &
            icore(iijka1), icore(iijka2),  &
            icore(iabij1), icore(iabij2),  &
            icore(itmp1), icore(itmp2), icore(itmp3), &
            icore(itmp4), icore(itmp5), icore(itmp6), &
            icore(iabci1), icore(iabci2), &
            icore(iend), e4t, e5st)
elseif(prec==8.and.totgpu==1) then
   write(6,*)'================calculate aab in gpu .dp================='
   if(iend > icrsizg) then
       write(6,*)'attention!!!! memory not enough for (T) calculations'
       write(6,*)'integrals would be divided into small part'
       igpu=1; igpu2=1
       call ccsdt_sym_aab_dg_l2(icore(ie), icore(it1aa), &
            icore(it21), icore(it22),  &
            icore(it221), icore(it222),   &
            icore(iijka1), icore(iijka2),  &
            icore(iabij1), icore(iabij2),  &
            icore(itmp1), icore(itmp2), icore(itmp3), &
            icore(itmp4), icore(itmp5), icore(itmp6), &
            icore(iabci1), icore(iabci2),  &
            icore(iend), e4t, e5st, icorg)
   else 
    call myicopy_g(icore(ie), icorg(ie), nva+noa+nvb+nob)
    call ccsdt_sym_aab_dg(icorg(ie), icorg(it1aa), &
            icorg(it21),  icorg(it22),    &
            icorg(it221), icorg(it222),   &
            icorg(iijka1),icorg(iijka2),  &
            icorg(iabij1),icorg(iabij2),  &
            icorg(itmp1), icorg(itmp2),  icorg(itmp3), &
            icorg(itmp4), icorg(itmp5),  icorg(itmp6), &
            icorg(iabci1),icorg(iabci2),  &
            icorg(iend), e4t, e5st, icore(it1aa))
  endif
elseif(prec==4.and.totgpu==1) then
   write(6,*)'================calculate aab in gpu .sp================'
   if(iend > icrsizg) then
       write(6,*)'attention!!!! memory not enough for (T) calculations'
       write(6,*)'integrals would be divided into small part'
       call ccsdt_sym_aab_sg_l2(icore(ie), icore(it1aa),  &
            icore(it21), icore(it22), &
            icore(it221), icore(it222),  &
            icore(iijka1), icore(iijka2),  &
            icore(iabij1), icore(iabij2),  &
            icore(itmp1), icore(itmp2), icore(itmp3), &
            icore(itmp4), icore(itmp5), icore(itmp6), &
            icore(iabci1), icore(iabci2),  &
            icore(iend), e4t, e5st, icorg)
   else
   igpu=1; igpu2=1
    call myicopy_g(icore(ie), icorg(ie), nva+noa+nvb+nob)
    call ccsdt_sym_aab_sg(icorg(ie), icorg(it1aa), &
            icorg(it21),  icorg(it22),    &
            icorg(it221), icorg(it222),   &
            icorg(iijka1),icorg(iijka2),  &
            icorg(iabij1),icorg(iabij2),  &
            icorg(itmp1), icorg(itmp2),  icorg(itmp3), &
            icorg(itmp4), icorg(itmp5),  icorg(itmp6), &
            icorg(iabci1),icorg(iabci2),  &
            icorg(iend), e4t, e5st, icore(it1aa))
  endif
endif
write(6,*)'e4t and e5st test for aab', e4ttest, e5sttest
write(6,'(A, F8.2, F8.2, F8.2)')'tgen,tvpq and ttmpe', tgen, tvpq, ttmpe
write(6,'(A, F8.2)')            'ttmpe2             ', ttmpe2
!
igpu=0; igpu2=0; totgpu = 0
call cpu_time(tt(3)); call walltimeso(rt(3))
1000 write(6,*)'==============finished CCSD(T) calculations=============='
write(6,*)'time for aaa', tt(2)-tt(1), rt(2)-rt(1)
write(6,*)'time for aab', tt(3)-tt(2), rt(3)-rt(2)
write(6,*)'time for all', tt(3)-tt(1), rt(3)-rt(1)
eccsdt = e4t*2+e5st*2
eccsdt = e4ttest*2+e5sttest*2
write(6,*)'E4T and E5ST      ', e4ttest*2, e5sttest*2
write(6,*)'CCSD(T) energy    ', eccsdt
if(ldebug) then 
   write(6,*)'E4T and E5ST test ', e4t*2, e5st*2
   write(6,*)'CCSD(T) test      ', e4t*2+e5st*2
endif
open(555, file='CCSD_ENERGY')
read(555,*) escf, eccsd
write(555,*) eccsdt
close(555)
write(6,*)'ESCF        ', escf
write(6,*)'ECCSD       ', eccsd
write(6,*)'ECCSD(T)    ', eccsdt
write(6,*)'Total ENERGY', escf+eccsd+eccsdt

!
!
return
end
